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PACS 32 . 80 . Qk - Coherent control of atomic interactions with photons 
PACS 32 . 80 . Fb - Photoionization of atoms and ions 

PACS 31.15.ac - High-precision calculations for few-electron (or few-body) atomic systems 

Abstract. - We demonstrate how the shape of femtosecond laser pulses can be tailored in order 
to obtain maximal ionization of atoms or molecules. For that purpose, we have overlayed a 
direct-optimization scheme on top of a fully unconstrained computation of the three-dimensional 
time-dependent Schrodinger equation. The procedure looks for pulses that maintain the same 
total length and integrated intensity or fluence as a given pulse that serves as an initial guess. 
It allows, however, for changes in frequencies - within a certain, predefined range - and overall 
shape, leading to enhanced ionization. We illustrate the scheme by calculating ionization yields 
for the molecule when irradiated with short (~ 5 fs), high- intensity laser pulses. 
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Introduction. — When atoms or molecules are irra- 
diated with laser fields that are intense enough to induce 
nonlinear effects, a wealth of fascinating phenomena may 
Oe observed [1]. This applies even to deceitfully "unin- 
teresting" systems such as the simplest molecule, H2" [2]: 
4bove-threshold or tunneling ionization [3], bond soften- 
ing [4] , bond hardening (light induced bound states or vi- 
brational trapping) [5] , charge resonance enhanced ioniza- 
tion [6,7], above threshold dissociation [8], high harmonic 
generation [9], etc. 

This very same complexity in the molecular reaction, 
however, is what permits to envision the possibility of con- 
trolling molecules with short (femtosecond time scale) and 
intense (10^^ — 10^^ W/cm^) laser pulses [10]. The short 
durations allow for coherent control: the systems evolve 
uncoupled to the environment, and can be steered towards 
the desired outcomes without relying in the more tradi- 
tional control parameters, i.e., average, thermodynamic 
functions such as the temperature. The high intensities 
trigger the strongly nonlinear, even non-pert urbative, re- 
sponse of the systems. An essential ingredient to realize 
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the molecular control is the capability of shaping the laser 
pulses - a technological area that has witnessed spectacu- 
lar advances in the recent years [11]. 

Yet this complexity implies the need for challenging the- 
oretical models. Not surprisingly, ionization is the first 
and most studied process, mainly because it could already 
be studied for atoms [12], and because in this intensity 
regime it almost always occurs, be it accompanied or not 
by other phenomena. Even in the absence of influence 
from nuclear dynamics, the ionization of molecules is sig- 
nificantly more complex than that of atoms, due to the 
electron emission from different atomic centers [13, 14]. 
Two rather successful models for molecular ionization that 
have recently been suggested are the so-called molecular- 
orbital strong-field approximation [13] and the molecular 
extension of the Ammosov-Delone-Kraine (ADK) approx- 
imation [15]. However, these approaches are insufficient 
as general tools [16]. 

A common feature of the approaches mentioned above 
is the use of the single-active-electron (SAE) approxima- 
tion. One-electron molecular systems such as the hydro- 
gen molecular ion H2" are therefore perfect candidates to 
isolate the error introduced by the SAE approximation 
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from further simphfications. Electron correlation orig- 
inates difficult and interesting phenomena such as non- 
sequential ionization [17]. 

In order to properly investigate the interaction of short 
and intense laser fields with molecules, one needs to per- 
form explicitly time-dependent calculations, even if it 
might imply a heavy computational burden. Calcula- 
tions of this kind, that propagate the time-dependent 
Schrodinger equation (TDSE), have been presented for H2" 
in the past, for example with the purpose of understanding 
the presence of maxima in the ionization yield for particu- 
lar internuclear separations [18], or in order to disentangle 
the relationship between ionization and dissociation [7,19]. 
Recently, Selst0 et al. [20] and Kjeldsen et al. [21] have re- 
ported calculations on the orientation dependence of the 
ionization yield - lifting the commonly used assumption 
of a molecular axis parallel to the light polarization. 

In this work we take a further step, and focus on the pos- 
sibility of theoretically designing, via fixed-nuclei three- 
dimensional (3D) TDSE calculations, laser pulses able to 
control (in particular, significantly enhance) the ioniza- 
tion yields, taking II2" as an example system. Some recent 
experimental breakthroughs on this area have triggered 
our interest. For example, Suzuki et al. [22] demonstrated 
the control of the multiphoton ionization channels of I2 
molecules by making use of a pulse shaping system capa- 
ble of varying in time the polarization directions. Simul- 
taneously, Brixner et al. [23] have made use of a similar 
polarization-shaping system to enhance ionization yields 
of diatomic molecules (K2). Our focus is, however, on lin- 
early polarized ultrashort pulses (~ 5 fs), so rapid that the 
nuclear movement does not play a role during the pulse 
action - in contrast to the studies in which the ionization 
is studied as the internuclear distance changes, leading to 
possible resonances. 

Methodology. — The optimization problem could be 
formulated in the language of quantum optimal-control 
theory (QOCT) [24]. It consists of a set of equations - 
along with various suitable, iterative algorithms that solve 
them - whose solution provides an optimized control field 
that typically maximizes a target operator O. In order to 
enhance ionization, one would just define O as the pro- 
jection onto unbound states, or, alternatively, the identity 
minus the projection onto the bound states: 



bound 



(1) 



However, we have experienced numerical difficulties when 
attempting to solve the QOCT equations for this partic- 
ular operator: The forward-backward propagations that 
must be performed in order to solve the QOCT equations 
proved to be, for our particular implementation, numeri- 
cally unfeasible when using the operator given in Eq. (1) 
to define the target. This was due to the appearance of 
fields with unrealistically high frequencies and/or ampli- 
tudes. We believe that the reason lies in the fact that 



the backward propagation must be performed after acting 
with the operator O on the previously propagated wave 
function. This eliminates the smooth, numerically friendly 
part of the wave function, enhancing, on the contrary, the 
high frequency components. This procedure is repeated at 
each iteration, eventually making the propagation impos- 
sible. We do not claim, however, that any other numerical 
implementation will not be able to successfully cope with 
this problem. 

Therefore, we have employed and present here, a direct 
optimization scheme, which is in fact much closer in spirit 
to the techniques utilized by the experimentalists [25]. In 
this scheme, we construct a merit function by considering 
the expectation value of the operator defined in Eq. (1) at 
the end of the propagation: 
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where x is the set of parameters that define the laser pulse, 
and l^a^(T)) is the wave function that results from per- 
forming the propagation with the laser determined by x, 
at the final time T. Of course, the sum over the bound 
states has to be truncated; for the calculations presented 
below, we find it sufficient to include the lowest ten states. 
The merit function is calculated by performing consecutive 
TDSE propagations: The resulting function values are fed 
into a recently developed derivative-free algorithm called 
NEWUOA [26] . This algorithm seeks the maximum of any 
merit function F{x) depending on N variables x, and does 
not necessitate the gradient VF. It is very effective for N 
larger than ten and smaller than a few hundred, which is 
the case considered here. In all our runs, we necessitated 
around two hundred iterations to converge the ionization 
yields to within 1%. 

The system and the TDSE propagations are modeled in 
our homegrown octopus code [27]. We represent the wave 
functions on a real-space rectangular regular grid, and fix 
the nuclear position at their equilibrium distance. The 
small length of the pulses used here justifies this simplifi- 
cation. We perform calculations setting the polarization 
direction both parallel to the molecular axis and perpen- 
dicular to it. The size of the simulation box is selected 
large enough to ensure that very little of the electronic 
density has reached the grid boundaries at the end of the 
laser pulse. Nevertheless, we add absorbing boundaries to 
remove this charge; if the propagation is pursued after the 
pulse, part of the density will "abandon" the simulation 
box; the remaining integrated density should approach (as 
it does) one minus the ionization probability calculated as 
the expectation value of Eq. (1). 

The laser pulse is taken in the dipolar approximation, 
and represented in the length gauge. The temporal shape 
of the pulse is given by a function /(t), which we expand 
in a Fourier series: 



fit) =fo + 



N 

E 

n=l 



fn\l 7^ COs(cJn^) ^ 9n\l S\Ii{uJnt) 



(3) 



Femtosecond laser pulse shaping for enhanced ionization 



with uon = 27rn/T. In order to ensure a physically mean- 
ingful laser pulse [28], we must have dtf{t) = 0, which 
implies /o = 0. Moreover, we must have /(O) = f{T) = 0, 
where T is the total propagation time. This poses the 
following constraint: 

N 

The sum over frequencies is truncated according to phys- 
ical considerations: Any pulse shaper must have a pre- 
defined range of frequencies it can work with. The feasi- 
bility of the numerical scheme depends on the possibility 
of truncating the previous expression at a reasonably low 
number. In the cases considered here, due to the short 
duration of the pulses, we obtain no more than around 20 
degrees of freedom by setting the maximum frequency to 
one Hartree. 

Evidently, by increasing the intensity of a pulse one can 
enhance the ionization yield. Our wish is to improve this 
yield by changing the pulse shape, and not simply by las- 
ing with larger intensity. Therefore, to ensure the fairness 
in the optimization search, we constrain the search to laser 
pulses whose time-integrated intensity (fluence) is prede- 
fined to some value Fq: 

n=l 

The search space {/n, ^n} is thus constrained to the hyper- 
sphere defined by the previous equation. But we must 
add the condition given by Eq. (4), which further restricts 
the search space to a hyper-ellipsoid. By performing the 
appropriate unitary transformation, this can be brought 
again into a hyper sphere: 

2N-1 

Fo=Y.^n- (6) 

This equal-fiuence condition can be guaranteed if we per- 
form a new transformation to hyper spherical coordinates: 

6 = Fo^/'cos(^i), 

6 = Fo'/'sin(^i)cos(^2), 

1/2 

^2N-2 = Fq' sin(6>i)...sin(6>2Ar-3)cos(6>2Ar-2), 

1/2 

i2N-i = Fq' sin(6>i)...sin(6>2Ar-3)sin(6>2Ar-2). (7) 

The set of angles {9j}^^^^ constitute the 2N — 2 vari- 
ables that define the search space for the optimization al- 
gorithm. 

Results. — The initial laser field before the optimiza- 
tion is a linearly polarized eight-cycle pulse having a si- 
nusoidal envelope, fixed peak intensity, and wavelength 
of A = 400 nm - a typical value for frequency-doubled 




peak intensity (10^^ W/cm^) 

Fig. 1: (color online). Ionization probability for the initial pulse 
(circles) and for the optimized pulse (squares) as a function of 
the peak intensity of the initial pulse. The polarization of the 
pulse is (a) parallel and (b) perpendicular to the molecular 
axis. 

Titanium-sapphire lasers. Correspondingly, the initial fre- 
quency is (jOq = 0.114 Ha, and the pulse length is 5.3 
fs. The maximum allowed frequency of the optimized 
pulse is set to cc^max = 2uoq. The pulse polarization is 
fixed to be parallel or perpendicular to the molecular 
axis. During the QOCT procedure, the polarization and 
the fluence Fq are kept fixed, but the peak intensity may 
change from the initial value, which is selected in the range 
/ = 0.5, 0.75, . . . , 2 X 10^^ W/cm^. 

Figure 1 shows the ionization probabilities as a function 
of the peak intensity (of the initial guess pulse) for the ini- 
tial and optimized pulses polarized parallel (a) and per- 
pendicular (b) to the molecular axis, respectively. Over- 
all, the pulse optimization leads to a significant increase in 
the ionization. As expected, the ionization yield is slightly 
larger for pulses polarized parallel to the molecular axis. 

To get more insight into the optimized ionization pro- 
cess, we plot in Fig. 2 the initial and optimal laser pulses 
and the occupations of some single-electron states during 
the pulse interaction. The peak intensity of the initial 
pulse is / = 2 X 10^^ W/cm^. Optimized pulses of both 
parallel (a) and perpendicular (c) polarization have large 
peaks near the end of the pulse. According to the corre- 
sponding occupations shown in Figs. 2(b) and (d), these 
amplitude peaks account for almost all of the ionization: 
During the peaks the ground-state occupations rapidly 
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Fig. 2: (color online), (a) Initial and optimized pulses (parallel 
polarization) and (b) the occupation of selected single-electron 
states in the optimized ionization process, when I — 2 x 10^^ 
W/cm^ (c-d) Same as (a-b) but for perpendicular polariza- 
tion. 



collapse. The 2p cj^ (2p tTu) excited state contributes to 
the process to a small extent in the parallel (perpendicu- 
lar) case, whereas the other states are involved by a nearly 
negligible fraction; overall, no excited bound states con- 
tribute significantly. Hence, within the constraints set here 
for the laser pulse, the optimal ionization of H2" is a direct 
process obtained by focusing most of the available pulse 
energy in a very short time frame - though keeping the 
integrated total field at zero in accordance with Maxwell's 
equations (see Ref. [28]). The electron densities during the 
ionization process are visualized in Fig. 3 for both parallel 
(a) and perpendicular (b) polarizations. 

The previous optimizations have produced rather "unin- 
teresting" solutions. It is a well known fact that in a short 
intense laser pulse, most of the ionization occurs during 
the peaks in the electric field. Therefore, the optimiza- 
tions have just attempted to create short, intense bursts 
of light. This fact can be understood if we consider the 
process happening in the quasi-static, tunneling regime, 
in which the total ionization can be approximated by con- 
sidering at each moment in time the static ionization rate 
that corresponds to the electric amplitude. This ioniza- 
tion rate is nonlinear, and it is much larger at the electric 
field peaks, which therefore cause most of the ionization. 
Note, however, that the cases discussed above lie in an 
intermediate regime between the tunneling and the multi- 
photon regime - the Keldysh parameter, 7, is of the order 
of one [The Keldysh parameter 7 is defined as yJ\Ei\/2Up^ 
where Ej is the ionization potential of the system, and Up 
is the pondemorotive energy, given in atomics units by 
(£^o/2cc;)^, ^0 being the peak intensity of the electric field, 
and uo the pulse frequency. Since our optimized lasers 
do not have a single frequency - not even necessarily a 
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Fig. 3: (color online), (a-f) Snapshots of the electron densities 
in the optimized ionization process when the pulse polarization 
is parallel to the molecular axis and the intensity is / = 2 x 10^^ 
W/cm^. The black dots mark the positions of the nuclei. The 
time-steps between snaphots are all equal, (g-1) The same for 
a pulse polarized perpendicular to the molecular axis. 

dominant one, we can only speak of approximate Keldysh 
parameters.] 

Alternatively, one can explain the simplicity of the 
pulses considering that the maximum allowed frequency, 
2(jOq = 0.228 Ha, is smaller than any resonance transition 
energy from the ground state. As a consequence, the sys- 
tem does not significantly populate these states, and the 
only ionizing channel is direct transition to the continuum. 

The picture changes significantly, however, if we allow 
for a larger cutoff frequency. First, this increases the 
value of the Keldysh parameter associated with the pro- 
cess, which may change the regime from a more quasistatic 
to a more of a multi-photon-like character. Secondly, the 
excited bound states are now accessible for single-photon 
transitions. For example. Fig. 4 displays results obtained 
for 4:UjQ. The intensity is here set to 0.5 x 10^^ W/cm^. 
Doubling the cutoff frequency of the search space has a 
significant effect in the total ionization yield: Now we ob- 
tain 0.99 for the ionization probability, whereas in the first 
optimization the yield was 0.20 (see Fig. 1, top panel, first 
point in the series). Note that the initial yield before any 
optimization was only 0.005. 

Moreover, the manner in which the ionization occurs 
with a larger cutoff frequency is qualitatively very differ- 
ent. Figure 4 displays the evolution of the occupation of 
some of the bound states. The first excited state (2p<Ju) 
plays a significant role, which can be understood because 
the transition energy from the ground state is now acces- 
sible in the field search space. In addition, a couple of 
other lowest states contribute in the ionization process in 
an ascending order as a function of time. It should be 
noted, however, that only the a orbitals, where the nodes 
are perpendicular to the polarization axis, participate in 
the transitions, and tt orbitals, for example, are not al- 
lowed due to a different symmetry. As a consequence of 
the involvement of several states in the ionization process, 
the structure of the optimized laser pulse shown in the 
upper panel of Fig. 4 is much more complicated than the 
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Fig. 4: (color online). Upper panel: Optimized laser pulse 
for the ionization when the cutoff frequency is 4a;o (see text) 
and the intensity is fixed to 0.5 x 10^^ W/cm^. Lower panel: 
Occupation of a few lowest states during the pulse interaction. 

single-burst fields obtained in the previous calculations. 

Conclusions. — In conclusion, we have shown 
with three-dimensional time-propagations of the time- 
dependent Schrodinger equation how the precise temporal 
shape of a short intense laser pulse may affect significantly 
the total ionization yield of H2" at fixed internuclear sep- 
arations. Moreover, we have employed a gradient-free op- 
timization technique to find the laser pulse that enhances 
ionization. This optimization can be constrained in differ- 
ent ways, accounting for the limitations of physical sources 
- not all frequencies and intensities are available, and not 
all possible shapes can be constructed with the state-of- 
the-art pulse shapers (although the technology improves 
at a phenomenal rate). The results will differ depend- 
ing on these constraints: The optimized laser pulse may 
be the single burst of electric field that one would ex- 
pect by considering a process in the tunnelling regime, 
or a field with a more complicated structure that drives 
the system through intermediate states - the ionization 
can be enhanced by resonant transitions. Whereas in the 
former case it would be easy to design intuitively pulses 
that maximize the ionization, in the latter an optimiza- 
tion algorithm such as the one presented in this work is 
necessary. 

* * * 
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